# !/usr/bin/env python
# -*- coding: utf-8 -*-
# loading packages
import numpy as np
import pandas as pd
from pandas import datetime

# data visualization and missing values
import matplotlib.pyplot as plt
import seaborn as sns  # advanced vizs
import missingno as msno  # missing values

# stats
from statsmodels.distributions.empirical_distribution import ECDF
from sklearn.metrics import mean_squared_error, r2_score

# machine learning
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import Lasso, LassoCV
from sklearn.model_selection import train_test_split, cross_val_score

def Regression():
    seed = 123
    data = pd.read_csv('E:\培训教程\python\机器学习数学基础\统计分析\回归分析\汽车价格预测\Auto-Data.csv', na_values='?')
    pd.set_option('display.max_rows', 500)
    pd.set_option('display.max_columns', 500)
    # print(data.columns)    # 查看数据的列
    # print(data.info())     # 查看数据的基本信息，数据类型等
    # print(data.describe()) # 数据更详细的内容，只能展示数值型的数据

    '''数据可视化模块，可以可视化数据，并观测缺失值分布'''
    #date_visible(data)

    '''
    通过可视化，发现数据中normalized-losses一列的缺失值较多    
    对缺失值的处理(对于缺失值少的列，可以直接删除；如果一列中缺失值较多，需要对缺失值进行合理填充)
    可以发现 80% 的 normalized losses 是低于200 并且绝大多数低于125
    一个基本的想法就是用中位数来进行填充，但是我们得来想一想，这个特征跟哪些因素可能有关呢？
    应该是保险的情况吧，所以我们可以分组来进行填充这样会更精确一些,使用相关特征的数据进行分组，        
    求缺失列数据的均值填充
    首先来看一下对于不同保险情况的统计指标
    '''
    #print(data.groupby('symboling')['normalized-losses'].describe())

    '''使用dropna()函数直接将缺失值很少的列所在的行删除掉'''
    data = data.dropna(subset=['num-of-doors', 'price', 'peak-rpm', 'horsepower', 'stroke', 'bore'])
    '''使用symboling分组之后对应的均值填充normalized-losses中对应的数据
    例如symboling分组为2的normalized-losses字段的均值为50，则在填充symboling为2对应的normalized-losses缺失的值的时候，使用50'''
    data['normalized-losses'] = data.groupby('symboling')['normalized-losses'].transform(lambda x: x.fillna(x.mean()))

    '''计算两两特征之间的关系'''
    #date_cormatrix(data)


    '''将长宽高三个变量使用体积一个特征代替；可以在实际中使用有替代性的一个特征代替多个特征'''
    data['volume'] = data.length * data.width * data.height
    '''
    增加体积特征后，删除长宽高三个特征，并且删除特征相关性特别高的一个变量，当两个特征的相关性特别高时，可以删除其中的一个
    保留一个即可，否则容易引起多重共线性(回归模型中两个或两个以上的自变量彼此相关的现象)
    分析数据之间的相关性之后，再原始数据集上直接操作，删除部分特征
    '''
    data.drop(['width', 'length', 'height', 'engine-size', 'curb-weight', 'city-mpg'], axis=1, inplace=True)
    #print(data)

    '''以上是使用pd.df对数据进行操作，观察特征之间的相关性之后对原始数据集的操作，也可以使用可视化图形观察特征之间的相关性'''
    '''绘制变量之间的相关关系图'''
    #heatmap(data)

    '''pairplot 绘制两两变量之间的散点图关系，
    hue='fuel-type' 按照fuel-type列的变量分类，palette设置调色板的颜色'''
    #pairplot(data)


    '''再来观察下车的价格和马力之间关系
    lmplot() 可以用来绘制线性线
    前三个参数分别为x轴，y轴，数据集，hue分类列，col 图形的列上的值，row图形的行的值
    此函数可以深入分析两个变量之间的线性关系'''
    sns.lmplot('price', 'horsepower', data, hue='fuel-type', col='fuel-type', row='num-of-doors', palette='plasma',
               fit_reg=True)
    #plt.show()

    features,target = data_sacler(data)
    '''划分数据集，将数据划分为训练集和测试集,其中测试数据战30%'''
    X_train, X_test, y_train, y_test = train_test_split(features, target, test_size=0.3, random_state=seed)
    # print("Train", X_train.shape, "and test", X_test.shape)

    '''Lasson回归
    # 相比线性归回，多添加了一个绝对值来惩罚过大的系数，alphas=0就和最小二乘法一样了'''
    alphas = 2. ** np.arange(2, 12)
    scores = np.empty_like(alphas)
    lassocv = LassoCV(cv=10, random_state=seed)  # 对数据进行交叉验证,cv表示需要切分的数据块数
    lassocv.fit(features, target)  # 使用特征数据预测价格
    lassocv_score = lassocv.score(features, target)
    lassocv_alpha = lassocv.alpha_  # 调用属性获取其alpha_的值
    print(lassocv_alpha)

    plt.figure(figsize=(10, 4))
    plt.plot(alphas, scores, '-ko')
    plt.axhline(lassocv_score)
    plt.xlabel('$alpha$')
    plt.ylabel('CV Score')
    plt.xscale('log', basex=2)
    sns.despine(offset=15)
    # It's already a good result: R squared of 0.82. Let's take a look at the features contributing to the model.
    print('CV results:', lassocv_score, lassocv_alpha)
    plt.show()

    model_l1 = LassoCV(alphas=alphas, cv=10, random_state=seed).fit(X_train, y_train)
    y_pred_l1 = model_l1.predict(X_test)
    print(model_l1.score(X_test, y_test))

    # 残差点
    plt.rcParams['figure.figsize'] = (6.0, 6.0)
    preds = pd.DataFrame({"preds": model_l1.predict(X_train), "true": y_train})
    preds["residuals"] = preds["true"] - preds["preds"]

    # 最好的情况下，残差点都在0附近浮动
    preds.plot(x="preds", y="residuals", kind="scatter", color=c)
    plt.show()
    print(MSE(y_test, y_pred_l1))
    print(R2(y_test, y_pred_l1))


def MSE(y_true, y_pred):
    mse = mean_squared_error(y_true, y_pred)
    print('MSE: %2.3f' % mse)
    return mse


def R2(y_true, y_pred):
    r2 = r2_score(y_true, y_pred)
    print('R2: %2.3f' % r2)
    return r2

def date_visible(data):
    '''
    对缺失数据的可视化展示，图形中黑色的表示有数据，白色的表示有缺失值
    '''
    msno.matrix(data)
    plt.show()
    '''从上面的图形中可以看出，normalized-losses 这一列缺失值却多,查看此列的数据'''
    print(data[pd.isnull(data['normalized-losses'])])

    # 将缺失值图形化，方便查看
    plt.figure(figsize=(12, 5))
    c = '#366DE8'
    # ECDF
    plt.subplot(121)
    cdf = ECDF(data['normalized-losses'])
    plt.plot(cdf.x, cdf.y, label="statmodels", color=c)
    plt.xlabel('normalized losses')
    plt.ylabel('ECDF')

    # 将缺失值删除之后，观察数据的分布情况，使用的是直方图
    plt.subplot(122)
    plt.hist(data['normalized-losses'].dropna(), bins=int(np.sqrt(len(data['normalized-losses']))), color=c)
    plt.show()

'''数据相关性'''
def date_cormatrix(data):
    '''特征相关性，会计算每个列与其他列的相关关系，返回的是一个以对角线对称的对称矩阵'''
    cormatrix = data.corr()
    #print(cormatrix)

    '''# 返回矩阵的上三角矩阵，并把对角线上的置位0，让他们不是最高的，目的是想取出来相关项最强的，
    # 按照相关性排序'''
    cormatrix *= np.tri(*cormatrix.values.shape, k=-1).T

    cormatrix = cormatrix.stack()
    print(cormatrix)
    print(type(cormatrix))

    '''按照Series(cormatrix)的数据的绝对值进行降序排序，reset_index()重新设置索引，原来的二重索引成为列'''
    cormatrix = cormatrix.reindex(cormatrix.abs().sort_values(ascending=False).index).reset_index()
    print(cormatrix)
    print(type(cormatrix))

    # 对DF的列进行重命令
    cormatrix.columns = ["FirstVariable", "SecondVariable", "Correlation"]
    print(cormatrix.head())

'''绘制两两变量之间的相关性'''
def heatmap(data):
    '''特征相关性，会计算每个列与其他列的相关关系，返回的是一个以对角线对称的对称矩阵'''
    corr_all = data.corr()
    mask = np.zeros_like(corr_all, dtype=np.bool)
    mask[np.triu_indices_from(mask)] = True
    f, ax = plt.subplots(figsize=(11, 9))

    '''绘制变量之间的热力图，表示相关性，颜色越深，表示相关性越强,当两个变量之间的相关性非常强的时候，可以删除其中的一个变量，只保留其中一个即可
    heatmap() 为热力图，不带参数也可以画出来
    参见： http://seaborn.pydata.org/generated/seaborn.heatmap.html?highlight=heatmap#seaborn.heatmap'''
    sns.heatmap(corr_all, mask=mask, square=True, linewidths=.5, ax=ax, cmap="BuPu")
    plt.show()

'''绘制两两变量之间的散点图'''
def pairplot(data):
    sns.pairplot(data, hue='fuel-type', palette='plasma')
    plt.show()

'''对数据进行标准化和归一化'''
def data_sacler(data):
    '''如果一个特性的方差比其他的要大得多，那么他可能支配目标函数，是估计不能像预期的那样正确的从其他特征中学习，
        这就是我们需要对数据进行缩放的原因，也成为数据标准化(以0为均值的正态分布)
        称之为数据的标准化和归一化'''

    '''对连续值进行标准化(以0为均值)'''
    target = data.price
    regressors = [x for x in data.columns if x not in ['price']]
    features = data.loc[:, regressors]   # 取出特定的列
    # 取出所有的数值型的列
    num = ['symboling', 'normalized-losses', 'volume', 'horsepower', 'wheel-base', 'bore', 'stroke','compression-ratio', 'peak-rpm']
    # scale the data 对数据进行缩放，不影响系数的评估
    standard_scaler = StandardScaler()
    features[num] = standard_scaler.fit_transform(features[num])  # 对数据进行标准化处理(基本上都要做)

    # 对分类数据进行one-hot编码，类似于对字符数据进行编码
    # 取出所有的字符类型变量
    classes = ['make', 'fuel-type', 'aspiration', 'num-of-doors', 'body-style', 'drive-wheels', 'engine-location','engine-type', 'num-of-cylinders', 'fuel-system']
    '''将分类变量转换为虚拟/指标变量'''
    dummies = pd.get_dummies(features[classes])
    '''# 将转换之后的列表拼接起来，并删除之前的特征，这回数据集中全都是数字了'''
    features = features.join(dummies).drop(classes, axis=1)
    return features,target

if __name__ == '__main__':
    Regression()